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SUMMARY 


This  technical  report  describes  a  computational  research  task  performed  to  support 
continuing  investigations  into  the  electromagnetic  behavior  of  photonic  crystals.  Simple 
dielectric  vein  structures  are  used  to  mimic  the  configuration  of  the  photonic  crystals  found 
within  the  wings  of  the  Green  Hairstreak  Butterfly,  an  insect  noted  for  its  use  of  photonic 
crystals  for  coloration.  An  earlier  investigation  employed  a  square  dielectric  vein  structure  to 
approximate  these  photonic  crystals.  In  this  case,  the  numerical  model  is  improved.  A 
rectangular  vein  structure  is  employed  with  two  distinct  values  for  the  vein  dielectric 
permittivity.  The  performance  of  the  rectangular  model  is  compared  with  that  of  the  square 
model  for  both  the  Transverse  Electric  and  Transverse  Magnetic  waves  under  mirror  symmetry. 
The  structure  of  attendant  photonic  bands  is  discussed. 
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1  INTRODUCTION 

1.1  Background 

The  application  of  optics  is  well  established  in  today’s  society  and  has  been  so  since  the 
days  of  Galileo’s  telescope  and  throughout  the  evolution  of  light  microscopy.  Optics,  the  physics 
of  light,  is  a  natural  subdivision  of  the  broader  physics  of  electrodynamics  and  involves  the  in- 
depth  study  of  electromagnetic  waves  propagating  through  various  media.  The  propagation  of 
electromagnetic  waves  can  be  controlled  by  manipulating  the  media  containing  the  light.  A  good 
example  of  this  idea  is  modern  day  optical  fiber.  An  optical  fiber  is  a  specially  designed 
waveguide  that  can  propagate  light  over  great  distances  with  little  or  no  losses.  This  type  of 
optical  device  and  attendant  waveguide  theory  is  discussed  more  thoroughly  in  Reference  [1], 

In  the  pages  that  follow,  we  turn  our  attention  to  a  different  type  of  optical  device  known 
as  the  Photonic  Crystal  (PC).  PCs  involve  a  higher  level  of  complex  light  wave  interactions  than 
do  optical  fibers.  In  one,  two  or  three  dimensional  space,  these  crystals  are  constructed  by  using 
a  periodic  distribution  of  dielectric  material.  A  dielectric  material  is  the  matter  carrying  the  light. 
Think  of  glass  or  even  water  as  examples  of  dielectric  media.  Like  optical  fibers,  PCs  propagate 
light  with  very  low  losses. [2]  The  periodic  distribution  of  dielectric  material  creates  an  array  of 
interfaces  in  the  crystal.  At  each  interface,  light  waves  both  scatter  and  transmit  electromagnetic 
energy.  The  superposition  of  incident  and  reflected  waves,  through  the  phenomenon  of 
interference,  may  “block”  the  propagation  of  light  of  certain  wavelengths.  If  light  of  all  possible 
propagation  directions  (wave  vectors)  for  a  finite  range  of  wavelengths  (or  frequencies)  is 
blocked,  the  crystal  exhibits  a  photonic  band  gap  (PBG).[2]  Band  gaps  can  exert  a  strong  effect 
over  the  colors  of  light  propagated  by  a  given  photonic  crystal.  These  crystals  are  purposeful  and 
commonly  exist  in  nature,  particularly  among  insects,  e.g.,  photonic  crystals  are  found  in  certain 
butterfly  wings  and  in  the  exo-skeletal  chitin  of  certain  species  of  beetle.  [3]  It  is  photonic 
crystals  that  are  believed  to  cause  the  brilliant  coloration  found  in  many  insects  and  birds.  [4] 

It  is  interesting  to  note  that  light  (photons)  interacting  with  the  dielectric  potential  exhibit 
many  of  the  same  characteristics  demonstrated  by  electrons  encountering  an  atomic  potential.  Of 
course,  the  atomic  potential  is  created  by  the  distribution  of  electrical  charge  existing  around 
atoms.  By  cleverly  tailoring  the  atomic  potential  within  a  tiny  silicon  structure,  one  may 
construct  a  semi-conductor  such  as  a  transistor.  This  type  of  device  is  governed  by  solid  state 
physics.  [5]  Because  of  the  similar  behavior  of  electromagnetic  waves,  we  can  apply  many  of  the 
same  techniques  used  in  solid  state  physics  to  study  photonic  crystals.  As  it  happens,  mainly  the 
differential  operator  and  boundary  conditions  change  when  migrating  from  the  study  of  electrons 
to  photons.  For  electrons,  we  utilize  the  Schrodinger  Equation  while  Maxwell’s  equations  are 
applied  for  photons.  [2] 

1.2  Objectives 

Section  2  begins  with  a  concise  exposition  of  the  practical  mathematics  surrounding  the 
solution  of  Maxwell’s  equations  for  light  waves  propagating  in  dielectric  media.  Our  analysis  is 
based  exclusively  in  the  frequency  domain  to  permit  the  identification  of  optical  frequencies  and 
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modes.  The  system  of  equations  is  placed  in  time -harmonic  form  to  derive  the  master  eigenvalue 
problem.  Vectors  within  the  crystal  lattice  are  then  used  to  derive  the  reciprocal  lattice  in  wave 
vector  space.  The  associated  Bloch  wave  vectors  are  used  to  place  the  master  equation  in  Bloch 
form. 


Section  3  contains  a  description  of  the  square  and  rectangular  crystal  lattices  of  interest. 
The  associated  reciprocal  lattices  are  presented,  and  the  key  points  outlining  the  Brillouin  zones 
for  each  reciprocal  lattice  are  shown.  The  distribution  of  dielectric  material  is  described  for  each 
2D  lattice. 

Section  4  presents  band  diagrams  (dispersion  relations)  for  the  square  and  rectangular 
lattices.  Two  different  values  of  the  relative  dielectric  permittivity  are  considered.  Also,  some  of 
the  mode  shapes  are  plotted  for  the  lower  frequency  Transverse  Electric  (TE)  and  Transverse 
Magnetic  (TM)  modes.  We  discuss  specific  differences  between  the  optical  modes  for  the  square 
and  rectangular  lattices.  Section  5  presents  a  summary  for  this  report  as  well  as  a  wrap-up  of  the 
conclusions. 
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2  THEORY 


2.1  Maxwell’s  Equations 

It  is  comforting  to  note,  but  perhaps  a  little  less  interesting,  that  even  though  many 
photonic  crystals  are  small,  we  may  treat  light  propagating  within  these  crystals  as 
electromagnetic  waves  instead  of  as  a  stream  of  photons.  That  is  to  say,  we  need  not  invoke 
quantum  optics  (yet). [6]  As  a  result,  this  report  will  regard  light  propagation  in  PCs  through  the 
use  of  classic  electromagnetism.  In  this  context,  the  physics  of  time-dependent  electromagnetic 
waves  is  governed  by  Maxwell’s  equations,  a  confluence  of  four  separate  partial  differential 
equations  derived  over  many  years  by  researchers  such  as  Coulomb,  Gauss,  Faraday  and 
Ampere. [7]  These  equations  were  finally  corrected  for  time  dependency  and  assembled  into  a 
system  by  J.C.  Maxwell. [8]  The  first  equation  in  this  system  is  Coulomb’s  (or  Gauss’)  Law,  i.e., 

V  *D=/?  (1) 

where  D ,  in  the  remainder  of  this  work,  is  denoted  as  the  displacement  field,  and  p  is  the  free 
charge  density,  a  function  defined  in  space.  By  assuming  that  space  contains  no  free  magnetic 
poles,  we  have  that  the  second  equation  is 


V  •  B  =  0  (2) 

(as  it  happens,  free  magnetic  monopoles  have  never  been  discovered),  where  B  is  denoted  the 
magnetic  induction  field.  The  third  law  is  due  to  Faraday,  i.e., 

V  xE  +  — -  =  0  (3) 

dt 


where  E  is  the  electric  field,  and  t  is  the  time  coordinate.  The  final  equation  is  Maxwell’s 
corrected  form  of  Ampere’s  Law.  One  may  recall,  that  in  its  original  form,  Ampere’s  “law”  was 
not  suitable  for  time  dependent  fields.  The  corrected  equation  is 


VxH-® 

dt 


(4) 


where  H  is  denoted  as  the  magnetic  field;  J  is  the  current  density  function.  These  equations  (1) 
through  (4)  are  remarkable  for  a  couple  of  different  reasons.  As  a  primary  consideration, 
Maxwell’s  equations  involve  vector-valued  quantities,  the  fields  B,  E,  H  and  D  as  well  as  the 
current  density  vector  J.  Secondly,  from  the  standpoint  of  mathematics,  we  have  eight  equations 
and  the  fifteen  unknowns  B,  E,  H,  D  and  J;  even  if  we  assume  that  p  is  known,  the  system  is 
not  closed.  Note  that  each  vector  component  is  a  scalar  unknown.  The  closure  difficulty  is 
resolved  by  invoking  constitutive  equations  to  relate  D  to  E  and  H  to  B.  We  have  that 


D  =  £-E 
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where  s  is  the  dielectric  permittivity  for  the  medium.  For  isotropic  media,  it  may  be  expressed 
in  terms  of  the  vacuum  permittivity  (the  permittivity  of  empty  space  ),  i.e.. 


or  s  =  s0  sr 


(6) 


Xe  is  the  medium’s  electric  susceptibility. [8]  The  term  sr  in  equation  (6)  is  denoted  as  the 

medium’s  relative  dielectric  permittivity.  In  some  cases,  this  material  property  is  frequency 
dependent.  For  anisotropic  media,  the  relationship  between  the  components  of  D  and  E  is  more 
complicated.  According  to  the  Bloembergen  model,  we  have  that 

=  £ijEj+XijkEjEk  (7) 


where  D;  and  Ej  represent  the  ith  components  of  the  displacement  and  electric  fields, 

respectively. [2]  Instead  of  the  simple  constants  used  in  (6),  siJ  and  %ijk  are  tensor-based 

quantities  that  represent  those  cases  where  the  point  dielectric  permittivity  changes  depending 
upon  the  direction  of  light  propagation. 

A  similar  constitutive  relationship  is  employed  to  relate  the  magnetic  and  magnetic 
induction  fields,  i.e., 


B  =  //H  (8) 

where  //  is  the  magnetic  permeability  of  the  light  propagating  medium.  This  parameter  is  often 
measured  with  respect  to  the  vacuum  permeability  ju0  .[8]  Using  a  linear  model  for  the 
permeability,  we  have 


A  =  Ao  0  +  lJ  (9) 

in  that  %m  *s  the  magnetic  susceptibility  of  the  medium.  [7]  For  anisotropic  magnetic  media,  we 
may  use  the  tensor-based  form 


Hi=EijBj  (10) 

a  slightly  different  form  of  (8). [8]  Normally,  the  magnetic  susceptibility  has  a  value  near  unity, 
but  it  may  take  on  other  values.  The  three  equations  represented  by  (10)  require  inversion  to 
match  the  mathematical  form  of  (8).  Equations  (5)  and  (8)  provide  six  additional  equations  to 
help  close  the  system.  A  potential  formulation  may  be  used  to  provide  a  final  closure  relation 
depending  upon  the  choice  of  solution  method  for  a  particular  problem. 

DISTRIBUTION  A 


4 


DISTRIBUTION  A 


2.2  Simplifying  Maxwell’s  Equations 

In  this  report,  our  goal  is  to  solve  Maxwell’s  equations  for  particular  two-dimensional 
photonic  crystals.  The  solution  is  performed  in  frequency  space;  both  light  frequencies  and 
optical  mode  shapes  are  to  be  determined.  In  order  to  enable  this  analysis,  Maxwell’s  equations 
are  first  cast  in  a  simplified  form.  We  begin  by  making  a  common  assumption  for  light  wave 
propagation  problems,  i.e.,  the  medium  contains  no  free  electric  charge  ( /?  =  0 )  and  no  current 
density  ( J  =  0 ).  The  dielectric  permittivity  is  assumed  to  be  independent  of  frequency  and  a 
function  of  position  only;  the  magnetic  permeability  is  set  at  the  vacuum  value.  The  space 
variables  for  the  field  are  annotated  by  the  position  vector  r.  With  these  simplifications,  the 
Maxwell  equations  can  be  written  as 


V«H(r,0=0 

(11) 

V.[s(r)E(r,0]  =  0 

(12) 

VxE(r,0  +  M)^r^  =  0 
ot 

(13) 

VxH(r,0^o^(r)^  =  0 

Ot 

(14) 

For  efficiency  of  notation,  the  relative  dielectric  permittivity  sr  is  now  and  hereafter  denoted  by 
s(v)  in  equations  (12)  and  (14). [2]  The  next  step  is  to  address  the  time  dependency  of  these 
equations. 

As  we  stated  earlier,  our  analysis  is  confined  to  frequency  space,  so  the  fields  are 
assumed  to  be  harmonic  (or  periodic)  in  time.  Hence,  they  may  be  expressed  in  the  form 

H(r,t)  =  H(r,t)exp(-icot) 

(15) 

E(r ,  t)  =  E(r,  t)  exp(-z  cot ) 

(16) 

where  co  is  the  angular  frequency  for  the  light,  and  of  course,  i 
Euler’s  formula,  exp(-icot)  =  cos(cot)-ism(cot).  When  (15)  and 
divergence  equations  (11)  and  (12),  we  obtain 

is  the  imaginary  unit.  Using 
(16)  are  substituted  into  the 

V  •  H(r)  =  0 

(17) 

V  •  [,c(r)  E(r)]  =  0 

(18) 
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The  Faraday  (13)  and  Ampere  (14)  equations  can  be  handled  in  a  similar  manner,  yet  an 
imaginary  factor  appears  when  the  complex  exponential  is  factored  from  these  equations.  The 
resulting  equations  are 

V  x  E(r)  -ia>ju0  H(r)  =  0  (19) 

VxH(r,l)  +  i®fof(r)E(r,O  =  0  (20) 

The  latter  two  equations  are  coupled,  but  by  following  Joannopoulos  [2],  they  may  be  decoupled. 
Equation  (20)  is  divided  by  the  relative  permittivity,  and  the  curl  of  the  result  is  taken  to  obtain 
the  intermediate  equation 


Vx 


'VxH(r)^ 
V  f(r)  y 


+  i  cos0  V  x  E(r)  =  0 


(21) 


Next  equation  (19)  is  solved  for  V  xE(r)  and  substituted  into  (21).  The  result  is 


Vx 


'VxH(r) 

V  £(r) 


\ 

) 


a2  s0  ju0  H(r)  =  0 


(22) 


Equation  (22),  denoted  as  the  master  equation,  is  quite  important;  note  that  c2  ,  the  speed  of  light 
squared  is  the  inverse  of  yjs0  ~ju0  ,  so  we  have  that 


Vx 


^VxH(r)^ 
v  £(r)  y 


f-1 

U  J 


H(r) 


(23) 


The  master  equation  consists  of  a  linear  differential  form  operating  on  the  magnetic  field  and  set 
equal  to  a  constant  times  the  magnetic  field  creating  an  eigenvalue  problem  with  eigenvalue 
(col cf  and  eigenvector  H(r).  The  linear  differential  form  (or  operator)  is  written  as 


0  =  V  x 


vA(r) 


\ 

Vx 

y 


(24) 


When  H(r)  has  been  determined,  the  electric  field  is  obtained  from  (20)  as 


E(r)  = - VxH(r)  (25) 

a>s0s(  r) 

2.3  Frequency  -  Wave  Vector  Solution 


The  primary  output  of  the  analyses  that  follow  is  a  numerical  relationship  between 
frequency  (co)  and  the  wave  vector  (k  =  (kx,ky,kz) )  for  propagating  light  waves.  This 
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relationship  is  often  referred  to  as  the  dispersion  relation.  The  wave  vector  magnitude  k  =  II  k  II  is 
simply  the  wavenumber  2n  I  A,  where  A  is  the  wavelength  for  the  light.  The  wave  vector 
implies  a  sense  of  direction  to  the  wave.  For  an  analysis  in  2D  or  3D,  there  is  a  locus  of  pertinent 
wave  numbers  revealed  by  the  theory.  For  example,  when  considering  a  2D  problem  at  the  plane 
of  mirror  symmetry,  we  can,  without  a  loss  of  generality,  assume  that  the  crystal  is  uniform 
along  the  z  axis  with  kv  =  0 .  As  a  result,  light  propagation  is  confined  to  the  xy-plane  with  the 
attendant  wave  vector  (k  =(kx,ky, 0)).  It  follows  that  the  dispersion  relation  is  a  locus  of  the 
form 


co  =  a>(kx,ky )  (26) 

where  kx  and  ky  vary,  in  the  most  general  sense,  between  0  and  qo.  More  specifically,  we  choose 

only  the  most  pertinent  wave  vectors  in  the  reciprocal  plane,  particularly  those  bordering  the 
Brillouin  zone. [2,9]  To  accomplish  this  task,  we  must  express  the  magnetic  field  in  Bloch  form. 

Bloch  form  is  closely  associated  with  the  idea  that  our  photonic  crystal  is  composed  of  an 
array  of  dielectric  material  that  is  periodic  (or  repetitive)  in  some  sense  across  the  entirety  of  the 
crystal  lattice.  More  specifically,  we  require  the  lattice  to  have  discrete  translational  symmetry 
with  a  lattice  constant  vector  a  .[2]  Let  the  vector  R  be  an  integer  multiple  of  the  lattice  constant 
vector,  i.e., 


R  =  na 


(27) 


Then  values  of  the  dielectric  permittivity  abide  by  the  relationship 

<?(r)  =  s(r  +  R)  (28) 

for  all  possible  integer  values  of  n .  Related  to  the  geometric  lattice  in  physical  space  is  the 
reciprocal  lattice  cast  in  wave  vector  space.  The  reciprocal  lattice  vectors  tu  are  defined  as 
follows. [10] 


b,  •  a,  =  In Stj  (29) 

Based  upon  this  definition,  the  Bloch  wave  vector  is  defined  in  3-space  as 

k  =  +  k2  b2  +  k3  b3  (30) 

Each  wave  vector  (30)  existing  within  the  Brillouin  zone  (of  the  reciprocal  lattice)  identifies  an 
eigenstate  with  angular  frequency  &>(k)  for  the  operator  (24). [2]  The  associated  eigenvector  (or 
eigenfunction)  has  the  form 


DISTRIBUTION  A 


7 


DISTRIBUTION  A 


Hk(r)  =  exp(/k  *r)uk(r) 


(31) 


In  this  eigenvector,  u*.  (r)  is  a  periodic  function  defined  on  the  physical  lattice.  We  can  find  a 

partial  differential  equation  for  this  function  by  substituting  (31)  into  the  master  equation  (23), 
i.e., 


Vx 


1 

(a>'\ 

V  xexp(ik  •r)uk(r) 

L*(r)  J 

U  J 

exp(/k  •r)uk  (r) 


d  fix. 

V (exp(/ k  •  r))  =  —  exp(i k  •  r )  =  exp(/ k  •  r )  i k .  — -  =  exp (i k  •  r)  ik.  on 


fix 


By  using  (33)  and  (34)  in  (32),  we  obtain 


Vx 


=  z'kexpO'k  *r) 


l 

^  co'' 

- exp(ik*r)  (Vxuk(r)  +  ikxuk(r)} 

_e(r) 

U  J 

exp(ik  •r)uk(r) 


Simplifying  (35)  further,  we  have  that 


Vx 


l 

r  co'' 

- exp(/k  *r)  {V  +  ik)xuk(r) 

l£(r)  ] 

U  J 

exp( ik  •  r)uk(r) 


(32) 


The  term  inside  of  the  square  brackets  can  be  simplified  by  applying  the  vector  identity 

Vx(/(r)g)  =  /(r)Vxg  +  V/(r)xg  (33) 

Equation  (33)  requires  V(exp(ik  «r)) .  This  expression  is  derived  as  follows. 


(34) 


(35) 


(36) 


Now  (33)  is  applied  again  to  the  square  bracketed  term  with  /(r)  =  exp(z‘k  »r)  and 
g  =  1/  £(r)(V  +  /k)uk(r) .  The  result  is 


exp(/k  *r)V  x 


(V  +  ik) 
e(r) 


xuk(r) 


+  Vexp(/k  «r)x 


(V  +  ik) 
s(r) 


xuk(r) 


(37) 


co 

\CJ 


exp(/k  *r)uk  (r) 

Continued  simplification  with  the  use  of  (34)  yields  a  new  eigenvalue  equation 
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/Y7  1  \  (V  +  /k)  1 

(V  +  ik)  x - — —  x  (uk(r)  = 

e(r) 


/  \2 

'  CO  A 


\C) 


uk(r) 


(38) 


The  operator  { •  }  in  (38)  is  also  Hermitian,  yet  it  is  not  used  alone. [2]  We  can  derive  a 
fundamental  constraint  equation  for  transverse  waves  by  substituting  the  Bloch  form  (31)  into 
(17).  With  the  use  of  the  vector  identity 


Vx(/(r)g)  =  /(r)V.g  +  V/(r).g 


(39) 


we  obtain  the  result 


( V  + )  k)  •  uk  (r)  =  0  (40) 

Equation  (38)  subject  to  the  constraint  (40)  and  the  periodicity  condition 

uk(r)  =  uk(r  +  R)  (41) 

can  be  solved  for  the  frequencies  and  eigenfunctions  describing  light  propagation  in  the  photonic 
crystal.[ll]  Solving  this  set  of  equations  on  a  dielectric  lattice  is  the  purpose  of  the  MIT 
Photonic  Bands  (MPB)  computer  program  applied  to  produce  the  results  shown  later  in  this 
report. 
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3  LATTICE  DESCRIPTIONS 

The  goal  of  the  present  work  is  to  present  and  contrast  optical  modes  for  a  set  of  square 
and  rectangular  dielectric  vein  structures.  These  2D  photonic  crystals  (uniform  in  the  third 
Cartesian  direction)  are  selected  for  study  because  of  their  likeness  to  a  simplified  view  of  the 
photonic  crystals  found  in  the  wings  of  Callophrys  rubi,  the  green  hairstreak  butterfly.  [3]  An 
electron  micrograph  of  a  part  of  the  wing  is  shown  in  Figure  1 .  While  the  computational  lattices 
shown  below  do  not  approach  the  level  of  complexity  of  the  actual  butterfly  wing,  some  of  the 
larger  geometric  structures  may  be  captured  by  an  elementary  numerical  model. 


Figure  1.  Electron  micrograph  of  photonic  crystals  located  in  the  wings  of  Callophyrs  rubi[ 3] 

3.1  Square  Dielectric  Vein  Lattice 

The  square  dielectric  vein  lattice  is  a  relatively  simple  2D  geometric  object.  [2]  The  lattice 
is  comprised  of  thin  horizontal  and  vertical  veins  of  dielectric  material  separating  square  regions 
filled  with  air.  Figure  2  contains  a  representative  picture  of  the  lattice  with  the  dielectric  veins 
(uniform  thickness)  colored  black;  air  filled  regions  are  white. 


Figure  2.  Diagram  of  the  Square  Dielectric  Vein  Lattice.  Unit  cell  is  repeated  three  times  in  the  horizontal  and  vertical 
directions 
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For  this  lattice,  the  physical  lattice  vectors  are: 


a  j  =  a  (1,0) 
a  2  =  a  (0,1) 


(42) 


Recall  that  Maxwell’s  equation  are  scale  invariant,  so  the  lattice  constant  “a”  in  (42)  is  an 
arbitrary  scaling  factor  (units  of  distance)  that  may  be  used  to  establish  a  lattice  of  any  desired 
size  as  long  as  we  do  not  violate  the  minimum  scale  for  the  electromagnetic  wave  model. 
Equation  (29)  may  be  used  to  show  that  the  reciprocal  lattice  vectors  are: 


2tt 

t>!  =  - (1,0) 

a 

2  n 

b2=—  (0,1) 
a 


(43) 


The  dielectric  veins  have  width  of  0.165a  (appropriately  scaled  with  the  lattice  constant).  Figure 
3  contains  a  plot  of  the  unit  cell  for  this  physical  lattice.  In  this  figure,  the  dielectric  material  is 


Figure  3.  Unit  cell  for  the  square  dielectric  vein  lattice.  Dielectric  boundaries  are  shown  in  red.  Unit  cell  boundary  is 
shown  in  green.  All  lengths  have  the  dimensions  of  the  lattice  constant  a 

outlined  with  the  color  red  and  has  a  cruciform  shape  centered  at  the  figure’s  center.  The 
boundary  for  the  unit  cell  is  delineated  by  the  green  line.  The  calculations  conducted  with  the  use 
of  this  lattice  geometry  are  dimensionless,  yet  the  x  and  y  coordinates  naturally  inherit  the 
dimensions  of  the  lattice  constant  a.  This  report  considers  two  values  for  the  relative  dielectric 
permittivity.  The  baseline  case  sets  s  =  8.9  to  provide  a  strong  jump  in  the  electric  field  at  the 
dielectric-air  interface.  The  baseline  is  selected  to  establish  a  dielectric/air  permittivity  difference 
of  nearly  10:1.  The  second  case  sets  s  =  1 .52  somewhat  closer  to  a  value  for  chitin,  a  material 
commonly  found  in  butterfly  wings.  In  future  studies,  s  is  to  be  raised  to  2.4,  a  common  value 
for  chitin. [12] 
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The  reciprocal  of  the  square  dielectric  vein  lattice  can  be  easily  produced  by  spanning  the 
reciprocal  lattice  vectors.  A  wave  vector  k  =  (kx,ky)  in  the  reciprocal  lattice  is  given  by  the 

formula 


k=/b1+mb2,  /,m  =  0, ±1, ±2, ... 


(44) 


More  specifically, 


k=  —  l,m  =  0,±1,±2,...  (45) 

a 

For  equation  (45),  dimensions  are  only  assigned  by  the  lattice  constant  a  ;  if  the  lattice  constant 
is  assigned  specific  length  dimensions,  k  then  possesses  dimensions  of  a  wave  number.  The 
reciprocal  lattice  for  the  square  dielectric  vein  configuration  is  shown  in  Figure  4. 


Figure  4.  Reciprocal  lattice  for  the  square  dielectric  vein  configuration.  Axes  are  delineated  with  dimensions  of  wave 
number.  The  red  points  shown  in  the  diagram  are  produced  by  the  mapping  (45).  The  red  square  outlines  the  Brillouin 
zone  for  this  lattice.  The  blue  lines  outline  the  irreducible  Brillouin  zone.  Key  points  are  defined  at  the  vertices  of  the 
irreducible  zone 

Figure  4  is  constructed  by  summing  integer  multiples  of  the  reciprocal  lattice  vectors.  The 
Brillouin  Zone  (BZ)  (more  appropriately  titled  as  the  first  Brillouin  zone)  is  produced  by  the 
procedures  described  in  References  2  and  9.  In  short,  it  is  minimum  area  figure  that  is  bounded 
by  parallel  rays  constructed  through  the  center  points  of  line  segments  extending  between 
reciprocal  lattice  points.  The  irreducible  Brillouin  zone  (IBZ),  shown  outlined  in  blue,  is  the 
minimum  area  locus  that  can  be  used  to  map  the  solution  in  the  rest  of  the  BZ  by  a  series  of 
simple  rotations  and  translations.  [2]  The  points  at  the  vertices  of  the  IBZ  are  key  locations  for 
calculation  of  the  frequency  space  solution. 

3.2  Rectangular  Dielectric  Vein  Lattice 

The  principal  photonic  crystal  of  interest  in  this  report  is  the  rectangular  dielectric  vein 
lattice.  The  square  lattice  discussed  in  the  preceding  section  is  actually  a  special  case  of  the 
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rectangular  lattice.  The  rectangular  lattice  is  of  interest  since  its  configuration  is  more  like  the 
array  of  chitin  shown  in  Figure  1 .  A  picture  of  the  rectangular  dielectric  vein  lattice  is  shown  in 
Figure  5.  Three  multiples  of  the  unit  cell  are  shown  in  each  direction. 


Figure  5.  Diagram  of  the  rectangular  dielectric  vein  lattice.  Unit  cell  is  repeated  three  times  in  the  horizontal  and  vertical 
directions.  Note  the  difference  between  vertical  and  horizontal  vein  widths 

In  Figure  5,  the  black  color  represents  dielectric  material  while  air  is  represented  by 
white.  For  this  photonic  crystal,  the  width  of  vertical  veins  is  0.13157a  ;  horizontal  veins  have  a 
width  of  0.052631a  where  a  is  the  lattice  constant.  The  air  gap  width  and  height  are  0.86842a 
and  0.21052a,  respectively.  These  parameters  are  chosen  to  match  the  aspect  ratio  of  the  chitin- 
air  cells  shown  in  Figure  1.  Two  parameters  are  required  to  describe  the  physical  lattice  vectors 
for  this  configuration.  These  vectors  are 


a1  =  a(1,0) 
a2  =  a  (0,0) 


In  this  case,  0  is  chosen  as  0.26315.  The  unit  cell  for  the  physical  lattice  is  shown  in  Figure  6. 
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Figure  6.  Unit  cell  for  the  rectangular  dielectric  vein  lattice.  Dielectric  boundaries  are  shown  in  red.  Unit  cell  boundary  is 
shown  in  green.  All  lengths  have  the  dimensions  of  the  lattice  constant  a 


The  reciprocal  lattice  vectors  are  derived  as 
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b1=— (1,0) 

a 

b2  =  —(0,1/0) 

a 


(47) 


Wave  vectors  in  the  reciprocal  lattice  are  defined  as  in  equation  (44)  with  the  specific  result 


k  = 


2n( , 


a 


l,—  ,  l,m  =  0,±1,±2,... 

V  V ) 


(48) 


A  diagram  of  the  reciprocal  lattice  for  the  rectangular  vein  structure  is  shown  in  Figure  7.  One 


kx 

Figure  7.  Reciprocal  lattice  for  the  rectangular  dielectric  vein  configuration.  Axes  are  delineated  with  dimensions  of  wave 
number.  The  magenta  points  shown  in  the  diagram  are  produced  by  the  mapping  (45).  The  blue  square  outlines  the 
Brillouin  zone  for  this  lattice.  The  red  lines  outline  the  irreducible  Brillouin  zone.  Key  points  are  defined  at  the  vertices  of 
the  irreducible  Brillouin  zone 


may  notice  that  the  IBZ  is  rectangular  in  shape  for  this  lattice.  [9]  For  the  square  lattice,  the  IBZ 
is  triangular.  A  possible  explanation  for  this  result  is  that  the  rectangular  lattice  lacks  a  symmetry 
occurring  in  the  square  lattice  (at  the  kx  =  ky  line).  This  result  is  very  interesting  and  begs  future 

investigation  concerning  symmetry  operations  in  the  Brillouin  zone.  The  lattice  geometries 
described  in  this  and  the  preceding  sections  summarize  the  input  used  for  calculating  the 
electromagnetic  behavior  of  the  square  and  rectangular  2D  photonic  crystals.  Again,  as  with  the 
square  lattice,  our  choice  of  relative  dielectric  permittivity  values  remains  the  same. 
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4  RESULTS 

The  electrical  field  results  (dispersion  relations  and  mode  shapes)  shown  below  have 
been  calculated  by  using  the  MIT  Photonic  Bands  (MPB)  computer  code  developed  by  the 
Joannopoulos  Ab  Initio  Physics  Group  at  the  Massachusetts  Institute  of  Technology.[ll]  The 
theory  behind  the  frequency  domain  solution  methods  applied  within  MPB  is  briefly  discussed  in 
Section  2  of  this  report.  In  this  section,  we  present  the  results  predicted  by  MPB  for  the  2D 
photonic  crystal  configurations  presented  in  the  preceding  section. 

4.1  Transverse  Electric  (TE)  Wave  Propagation  at  s  =  8.90 

Electromagnetic  wave  solutions  for  the  square  and  rectangular  vein  dielectric  lattices 
have  been  calculated  at  a  relative  dielectric  permittivity  ( s )  value  of  8.90.  Figure  8  contains  the 
dispersion  relations  (band  diagrams)  for  TE  wave  propagation  at  s  =  8.9  for  both  lattices. 


Figure  8.  Dispersion  relations  for  transverse  electric  (TE)  wave  propagation  in  the  square  and  rectangular  dielectric  vein 
lattices  at  relative  permittivity  8.9.  The  dispersion  relation  for  the  square  lattice  is  shown  in  the  left  plot  while  the  plot  for 
rectangular  lattice  is  presented  on  the  right.  The  notation  wave  vector  number  is  used  to  indicate  the  location  of  the  wave 
vector  on  the  boundary  of  the  irreducible  Brillouin  zone.  Only  the  first  five  eigenmodes  are  displayed 

Within  the  square  lattice,  a  band  gap  is  evident  for  TE  waves.  The  rectangular  lattice  exhibits  no 
band  gaps.  It  is  interesting  to  compare  2D  contour  plots  of  the  eigenmode  shapes  for  these  cases. 
We  begin  by  examining  the  TE  mode  shapes  for  s  =  8.9 ,  a  baseline  case  useful  for  comparison. 
Recall  that  for  the  TE  mode,  the  electric  field  is  confined  to  the  xy-plane;  its  z  component  is  zero. 
The  associated  magnetic  field  is  aligned  with  the  z  axis.  All  of  the  mode  shape  plots  that  follow 
are  computed  at  wave  vectors  (k-points)  k  =  2nla  (0.5,  0.5,  0)  for  the  square  lattice  and  at  k  = 
2n, la  (0.5,  1.9,  0)  for  the  rectangular  lattice.  These  k-points  are  identified  in  Figure  8  with  wave 
vector  number  1 1  for  both  lattices.  For  the  fundamental  frequency,  Figure  9  contains  the  plots  of 
the  real  part  of  Ex  for  each  lattice.  Note  that  both  lattices  are  graphed  on  plots  of  the  same  size 

even  though  their  aspect  ratios  differ.  Doing  so  facilitates  a  side-by-side  comparison  of  the 
electric  field  components.  As  is  evident  in  Figure  9,  the  difference  in  this  field  component  is 
quite  large  when  moving  between  the  square  lattice  (left)  and  the  rectangular  lattice  (right).  In 
fact,  Ex  is  many  times  stronger  in  the  rectangular  lattice  than  it  is  in  the  square  lattice.  The  mode 
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Figure  9.  Real  part  of  Ex  for  the  fundamental  frequency  computed  at  (left)  k-point  (wave  vector)  2 n/a  (0.5, 0.5,0)  in  the 
square  lattice  and  at  (right)  k-point  2 nla  (0.5, 1.9,0)  in  the  rectangular  lattice.  Vein  dielectric  permittivity  8.90 
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Figure  10.  Real  part  of  Ey  for  the  fundamental  frequency  computed  at  (left)  k-point  (wave  vector)  2 nla  (0.5, 0.5,0)  in  the 
square  lattice  and  at  (right)  k-point  2 nla  (0.5, 1.9,0)  in  the  rectangular  lattice.  Vein  dielectric  permittivity  8.90 


shape  is  more  widespread  on  the  rectangular  lattice  than  on  the  square  lattice.  On  the  square 
lattice,  variations  in  Ex  are  mostly  confined  to  the  vein  junctions  (or  spots  [2]).  At  these 

junctions,  Ex  exhibits  quadrupolar  variation.  For  the  rectangular  lattice,  Ex  is  significantly 
compressed  elongating  it  in  the  jc  direction.  It  is  also  interesting  to  see  that  nodes  are  located 
along  the  horizontal  veins.  Also,  the  periodic  variation  in  Ex  extends  well  beyond  the  lattice 
shown.  In  fact,  one  needs  to  extend  the  lattice  over  eight  unit  cells  in  order  to  see  the  field 
variation.  The  real  part  of  the  y  component  of  the  electric  field  E  is  shown  in  Figure  10.  In  this 
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case,  field  component  oscillations  for  the  square  lattice  are  still  centered  around  the  vein 
junctions,  but  possess  a  dipolar  character.  Nodes  are  situated  within  the  air  space  throughout  the 
lattice.  For  the  rectangular  lattice,  the  field  has  greater  extent  with  periodicity  extending  beyond 
the  section  of  the  lattice  shown.  The  extent  of  the  field  can  be  observed  if  more  cells  are  plotted. 
The  locus  of  field  nodes  occurs  along  horizontal  lines  through  the  center  of  the  air  space.  The 
peaks  in  Ey  are  also  out  of  phase  with  respect  to  the  peaks  in  Ex. 

We  continue  our  examination  of  TE  mode  propagation  by  considering  the  next  higher 
frequency  for  both  lattices.  Figure  1 1  contains  plots  of  the  real  part  of  Ex  at  the  same  k-points 
selected  for  the  square  and  rectangular  lattices.  It  is  evident  that  the  field  behavior  is  different  in 


SQUARE  RECTANGULAR 

Figure  11.  Real  part  of  Ex  for  the  second  frequency  computed  at  (left)  k-point  (wave  vector)  2 nla  (0.5, 0.5,0)  in  the  square 
lattice  and  at  (right)  k-point  In! a  (0.5, 1.9,0)  in  the  rectangular  lattice.  Vein  dielectric  permittivity  8.90 
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Figure  12.  Real  part  of  Ey  for  the  second  frequency  computed  at  (left)  k-point  (wave  vector)  2 nla  (0.5, 0.5,0)  in  the  square 
lattice  and  at  (right)  k-point  In! a  (0.5, 1.9,0)  in  the  rectangular  lattice.  Vein  dielectric  permittivity  8.90 
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this  case,  especially  for  the  square  lattice.  Figure  12  contains  the  associated  contour  plots  for  Ey. 
For  the  square  lattice  at  this  the  next  to  lowest  frequency,  the  differences  in  electric  field 
behavior  are  very  pronounced.  For  the  fundamental  frequency,  the  field  activity  is  concentrated 
near  the  vein  lattice  spots  (vein  junctions)  with  some  dipolar  activity  on  the  vertical  veins.  Since 
the  veins  possess  high  relative  dielectric  permittivity  £•  =  8.9,  we  expect  that  the  field  will 
concentrate  near  the  veins. [2]  For  the  second  frequency  we  see  different  behavior  for  the  square 
lattice.  Ex  nows  concentrates  in  the  air  space  between  the  veins,  i.e.,  in  the  air  band.  The  field 
variation  also  indicates  dipolar  variation  (horizontal  for  Ex  and  vertical  for  Ey).  Why  does  the 
electric  field  change  so  dramatically?  The  answer  is  simple.  Figure  8  clearly  shows  that  a  band 
gap  exists  on  the  square  lattice  for  the  TE  mode  between  the  fundamental  and  second 
frequencies. 

Differences  calculated  on  the  rectangular  lattice  when  transitioning  between  the  two 
frequencies  are  not  as  pronounced.  On  this  lattice,  Ex  still  has  nodes  located  along  horizontal 
lines  in  the  air  space.  The  highest  magnitude  values  of  Ex  lie  on  the  horizontal  veins.  Ey  is  the 
stronger  component  of  the  two.  This  behavior  may  be  caused  by  the  fact  that  horizontal  dielectric 
veins  are  closer  together  and  tend  to  concentrate  field  lines  nearby.  On  the  other  hand,  the  mode 
shape  (or  field  distribution)  for  Ey  changes  to  a  greater  degree.  The  peaks  and  troughs  in  this 
field  component  clearly  occur  in  the  air  spaces,  and  the  component  varies  more  slowly  in  the 
horizontal  direction  than  in  the  vertical  direction.  The  phase  difference  between  Ex  and  Ey 
remains  and  possesses  roughly  the  same  magnitude  as  for  the  fundamental  mode.  Of  course,  the 
second  frequency  is  very  close  to  the  fundamental  for  the  wave  vector  in  question  (at  wave 
vector  number  1 1  in  Figure  8).  Hence,  the  modes  behave  in  a  similar  manner. 


Continuing  with  our  investigation,  we  turn  our  attention  to  the  third  frequency  exhibited 
at  the  k-points  in  question.  Contours  of  the  real  part  of  the  x  component  of  the  electric  field  are 
shown  in  Figure  13  for  both  the  square  and  rectangular  lattices.  For  this  case,  field  variation  is 
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Figure  13.  Real  part  of  Ex  for  the  third  frequency  computed  at  (left)  k-point  (wave  vector)  2n/a  (0.5, 0.5,0)  in  the  square 
lattice  and  at  (right)  k-point  In! a  (0.5, 1.9,0)  in  the  rectangular  lattice.  Vein  dielectric  permittivity  8.90 
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again  centered  at  spots  in  the  square  lattice.  For  the  rectangular  lattice,  variation  in  51(£'x)is 
similar  to  that  exhibited  at  the  second  frequency.  Extrema  exist  at  loci  along  horizontal  veins  and 
at  the  centers  of  the  air  spaces.  Figure  14  contains  the  corresponding  plots  for  '.KtF  j .  On  the 


Figure  14.  Real  part  of  Ey  for  the  third  frequency  computed  at  (left)  k-point  (wave  vector)  2 nla  (0.5, 0.5,0)  in  the  square 
lattice  and  at  (right)  k-point  In! a  (0.5, 1.9,0)  in  the  rectangular  lattice.  Vein  dielectric  permittivity  8.90 


square  lattice,  the  variation  of  Ey  is  similar  to  that  of  Ex\  these  field  oscillations  are  concentrated 
again  around  the  vein  junctions  (spots)  and  possess  some  of  the  characteristics  of  quadrupoles. 
The  variation  of  Ex  at  the  same  frequency  is  more  dipolar.  On  the  rectangular  lattice,  the  field 
oscillations  are  still  predominantly  oriented  in  the  horizontal  direction.  The  field  nodes  exist 
horizontally  along  veins  and  at  the  air  gap  centers.  The  periodic  nature  of  the  field  extends  well 
beyond  the  section  of  the  lattice  shown,  and  the  phase  difference  is  still  exhibited  between  x  and 
y  electric  field  components. 

4.2  Transverse  Electric  (TE)  Wave  Propagation  at  e  =  1 .52 

In  this  section  of  the  report,  we  repeat  the  analyses  perform  in  Section  4.1  now  for  the 
relative  permittivity  of  £-  =  1.52.  The  dispersion  relations  for  the  square  and  rectangular 
dielectric  vein  lattices  are  shown  in  Figure  15.  The  first  characteristic  of  the  band  diagrams 
shown  is  that  neither  lattice  possesses  a  bad  gap.  This  result  is  anticipated  since  the  difference 
between  the  relative  permittivity  (index  of  refraction  squared)  for  the  dielectric  and  air  is  so 
small.  Secondly,  the  band  diagram  for  the  rectangular  lattice  is  sharply  peaked  and  jagged  at 
several  wave  vector  numbers,  far  more  so  than  for  the  square  lattice.  This  jaggedness  is  very 
pronounced  at  the  lower  frequencies  for  a  wide  range  of  wave  vectors.  The  rectangular  lattice 
also  reaches,  for  a  given  eigenmode,  a  higher  range  of  frequencies  than  does  the  square  lattice. 
These  regions  in  the  dispersion  relation  graph  correspond  to  the  edges  of  the  Brillouin  zone 
where  ky  changes  with  kx  is  fixed.  In  the  crystal  lattice,  this  region  corresponds  to  the  edges  that 
change  in  the  y  direction.  The  higher  frequencies  encountered  in  this  region  may  be  explained  as 
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follows.  In  order  to  fit  a  fixed  number  of  wavelengths  within  this  distance,  a  comparatively 
shorter  wavelength  and  thus  a  higher  frequency  must  be  used  (when  compared  to  the  square 


Figure  15.  Dispersion  relations  for  transverse  electric  (TE)  wave  propagation  in  the  square  and  rectangular  dielectric  vein 
lattices  at  relative  permittivity  1.52.  The  dispersion  relation  for  the  square  lattice  is  shown  in  the  left  plot  while  the  plot 
for  rectangular  lattice  is  presented  on  the  right.  We  notation  wave  vector  number  is  used  to  indicate  the  location  of  the 
wave  vector  on  the  boundary  of  the  irreducible  Brillouin  zone.  Only  the  first  five  eigenmodes  are  displayed. 


lattice).  When  examining  the  mode  shapes  for  this  case,  we  use  the  same  k-points  specified  in  the 
preceding  section. 

Figure  16  contains  TE  mode  shapes  for  the  x  component  of  the  electric  field  Ex  at  the 
fundamental  frequency  calculated  at  wave  vector  number  11  for  both  the  square  (left)  and 
rectangular  (right)  lattices.  For  the  square  lattice,  Ex  clearly  exhibits  dipolar  behavior  with  the 
field  nodes  aligned  along  the  vertical  veins.  Adjacent  horizontal  lattice  tiers  alternate  in  phase. 


Figure  16.  Real  part  of  Ex  for  the  fundamental  frequency  computed  at  (left)  k-point  (wave  vector)  2nla  (0.5, 0.5,0)  in  the 
square  lattice  and  at  (right)  k-point  2 nla  (0.5, 1.9,0)  in  the  rectangular  lattice.  Vein  dielectric  permittivity  1.52 
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On  the  rectangular  lattice,  Ex  shares  the  same  dipolar  character,  but  the  field  oscillations  are  more 
extensive  on  the  lattice.  The  nodes  are  oriented  along  vertical  aligned  loci  located  at  the  centers 
of  air  gaps.  These  vertical  loci  are  spaced  apart  at  a  periodicity  (in  terms  of  unit  cells)  exceeding 
the  extent  lattice  segment  shown  in  the  figure.  The  magnitude  of  the  variations  in  Ex  for  the 
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Figure  17.  Real  part  of  Ey  for  the  fundamental  frequency  computed  at  (left)  k-point  (wave  vector)  2 nla  (0.5, 0.5,0)  in  the 
square  lattice  and  at  (right)  k-point  2 nla  (0.5, 1.9,0)  in  the  rectangular  lattice.  Vein  dielectric  permittivity  1.52 
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Figure  18.  Real  part  of  Ex  for  the  second  frequency  computed  at  (left)  k-point  (wave  vector)  2 nla  (0.5, 0.5,0)  in  the  square 
lattice  and  at  (right)  k-point  (0.5, 1.9,0)  in  the  rectangular  lattice.  Vein  dielectric  permittivity  1.52 

rectangular  lattice  are  considerably  larger  than  those  calculated  for  the  square  lattice.  The  peak 
values  in  Ex  are  located  at  air  gap  centers.  Figures  17  present  the  y  component  of  the  electric  field 
for  the  fundamental  frequency.  This  component’s  behavior  is  interesting  because  on  the  square 
lattice,  it  possesses  almost  no  detectable  magnitude.  Hence,  the  electric  field  for  the  square 
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lattice  is  fully  aligned  with  the  x  axis.  On  the  other  hand,  the  rectangular  lattice  Ey  has  its  peak 
magnitudes  centered  along  the  horizontal  veins  with  nodes  situated  along  horizontal  loci  located 
at  the  air  gap  centers.  The  phase  difference  for  the  x  and  y  components  persists  as  shown  in  the 
previous  cases.  The  structure  of  the  electric  field  for  the  rectangular  lattice  exhibits  a  great  deal 
of  commonality  for  different  frequencies.  Consider  the  behavior  of  the  electric  field  for  the  two 
lattices  at  the  second  modal  frequency.  Ex  for  the  square  lattice,  shown  in  Figure  18,  exhibits 
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Figure  19.  Real  part  of  Ey  for  the  second  modal  frequency  computed  at  (left)  k-point  (wave  vector)  2n/a  (0.5, 0.5,0)  in  the 
square  lattice  and  at  (right)  k-point  2 nla  (0.5, 1.9,0)  in  the  rectangular  lattice.  Vein  dielectric  permittivity  1.52 

fluctuations  centered  at  the  vein  junctions  (spots)  and  within  the  air  gaps.  A  quadrupole-like 
variation  is  observed  within  the  air  gaps  while  between  vein  spots,  a  dipolar  fluctuations  exists. 
Ex  for  the  rectangular  lattice  exhibits  the  same  dipolar  variation  with  the  field  peaks  located 
along  the  horizontal  veins.  The  horizontal  node  loci  are  found  at  the  air  gap  centers.  Figure  19 
contains  the  plots  of  Ey  for  the  second  modal  frequency.  The  square  lattice’s  field  is  dominated 
by  dipolar  variations  that  recur  throughout  the  crystal.  These  variations  are  centered  on 
horizontal  veins  with  an  additional  lobe  of  the  fluctuation  located  at  the  air  gap  center.  For  the 
rectangular  lattice,  Ey  still  possesses  a  dipolar  shape  except  peak  field  magnitudes  are  located 
along  horizontal  loci  located  at  air  gap  centers.  The  field  contours  do  exhibit  some  higher  order 
characteristics  like  multiple  local  maxima  and  minima.  We  continue  with  our  investigation  by 
considering  electric  field  behavior  at  the  third  modal  frequency. 

Figure  20  contains  plots  of  Ex  for  the  third  modal  frequency  computed  for  the  square  and 
rectangular  lattice  structures  ( ).  The  square  lattice  is  characterized  by  the  presence  of  monopoles 
centered  at  the  vein  spots.  The  rectangular  lattice  exhibits  a  more  complicated  system  of  dipoles 
aligned  with  the  horizontal  veins.  The  peaks  are  centered  on  both  the  horizontal  veins  and  at  the 
air  gap  centers.  The  maximum  positive  and  negative  peaks  alternate  from  vein  to  air  gap  when 
scanning  in  the  y  direction  along  the  lattice.  The  plots  for  Ey  are  provide  in  Figure  21.  In  this 
case,  the  square  lattice  is  characterized  by  the  presence  of  more  monopoles,  this  time  located  at 
the  air  gap  centers.  The  positive  and  negative  peak  magnitudes  alternate  between  adjacent 
monopoles.  In  this  case,  nodes  form  loci  along  the  veins.  The  rectangular  lattice  clearly 
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Figure  20.  Real  part  of  Ex  for  the  third  modal  frequency  computed  at  (left)  k-point  (wave  vector)  2n/a  (0.5, 0.5,0)  in  the 
square  lattice  and  at  (right)  k-point  Inla  (0.5, 1.9,0)  in  the  rectangular  lattice.  Vein  dielectric  permittivity  1.52 


Figure  21.  Real  part  of  Ey  for  the  third  modal  frequency  computed  at  (left)  k-point  (wave  vector)  Inla  (0.5, 0.5,0)  in  the 
square  lattice  and  at  (right)  k-point  2 nla  (0.5, 1.9,0)  in  the  rectangular  lattice.  Vein  dielectric  permittivity  1.52 


exhibits  more  complex  dipole  behavior  with  peak  values  existing  along  horizontal  veins  and  at 
air  gap  centers.  In  this  case,  a  cycle  requires  2/3  of  a  unit  cell  (vertically)  for  completion.  As  in 
the  previous  cases,  there  is  an  outstanding  phase  difference  between  Ex  and  Ey.  It  seems  evident 
that  the  field  behavior  on  the  rectangular  lattice  is  governed,  to  a  great  degree,  by  field 
compression  created  by  the  reduced  spacing  between  horizontal  veins.  The  spacing  between 
horizontal  veins  is  roughly  26%  of  the  spacing  between  vertical  veins. 
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4.3  Transverse  Magnetic  (TM)  Wave  Propagation  at  s  =  8.90 


Our  discussion  of  electromagnetic  wave  propagation  on  the  two  test  lattices  is  incomplete 
without  including  a  discussion  of  TM  wave  propagation.  Recall  that  since  the  magnetic  field 
components  are  confined  to  the  xy-plane,  there  is  only  one  electric  field  component,  Ez ,  aligned 
in  the  z  direction.  Band  diagrams  have  been  prepared  for  TM  propagation.  Figure  22  contains 
these  diagrams  for  the  test  lattices  at  s  =  8.90 .  The  square  lattice  possesses  a  small  band  gap  of 
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Figure  22.  Dispersion  relations  for  transverse  magnetic  (TM)  wave  propagation  in  the  square  and  rectangular  dielectric 
vein  lattices  at  relative  permittivity  8.90.  The  dispersion  relation  for  the  square  lattice  is  shown  in  the  left  plot  while  the 
plot  for  rectangular  lattice  is  presented  on  the  right.  The  notation  wave  vector  number  is  used  to  indicate  the  location  of 
the  wave  vector  on  the  boundary  of  the  irreducible  Brillouin  zone.  Only  the  first  five  eigenmodes  are  displayed. 


width  5.3%  between  bands  3  and  4.  This  gap  is  indicated  by  the  yellow  bar  shown  in  Figure  22. 
The  MPB  computer  code  predicts  two  other  band  gaps,  but  these  gaps  are  so  small,  their 
existence  is  questionable,  i.e.,  they  spuriously  result  because  of  numerics.  The  rectangular  lattice 
contains  no  band  gaps,  but  for  larger  dielectric  permittivity  values,  one  may  form  between  modal 
frequencies  four  and  five.  The  rectangular  plot  is  characteristically  jagged  as  are  its  predecessors 
for  the  TE  propagation  mode.  We  can  now  present  the  TM  mode  shapes. 


At  the  fundamental  frequencies,  mode  shapes  for  Ez  are  compared  in  Figure  23  for  the 
two  lattices.  The  square  lattice  is  dominated  by  an  array  of  monopoles,  alternating  in  sense  and 
centered  at  the  vein  spots.  The  rectangular  lattice  possesses  dipole  characteristics  with  peak 
values  lying  along  horizontal  veins.  The  nodes  fall  in  two  regions.  The  first  is  located  along 
horizontal  lines  of  air  gap  centers.  The  second  region  lies  along  vertical  lines  of  gap  centers;  the 
periodicity  for  these  loci  extend  outside  of  the  section  of  the  lattice  shown.  The  peak  magnitudes 
of  Ez  differ  by  less  than  a  factor  of  three.  Figure  24  contains  plots  of  the  real  part  of  E-  for  the 
second  modal  frequency.  In  this  case,  the  square  lattice  is  characterized  by  four-lobed  wave 
structures  centered  at  the  vein  spots.  The  polarity  of  these  structures  alternates  between  adjacent 
spots.  Nodes  exist  at  the  air  gap  centers  and  at  vein  midpoints.  The  magnitude  of  these 
fluctuations  is  quite  small.  The  rectangular  lattice  retains  dipolar  behavior  with  a  horizontal 
periodicity  that  extends  beyond  the  lattice  segment  shown.  The  peak  field  magnitudes  for  the 
rectangular  lattice  exceed  those  for  the  corresponding  square  lattice  mode.  In  the  y  direction, 
peaks  in  E-  alternate  between  air  gap  centers  making  for  distinct  y-oriented  dipolar  wave 
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Figure  23.  Real  part  of  Ez  for  the  fundamental  TM  mode  computed  at  (left)  k-point  (wave  vector)  2n/a  (0.5, 0.5,0)  in  the 
square  lattice  and  at  (right)  k-point  2 nla  (0.5, 1.9,0)  in  the  rectangular  lattice.  Vein  dielectric  permittivity  8.90 


Figure  24.  Real  part  of  Ez  for  the  TM  mode  at  the  second  modal  frequency  computed  at  (left)  k-point  (wave  vector)  2 nla 
(0.5, 0.5,0)  in  the  square  lattice  and  at  (right)  k-point  2nla  (0.5, 1.9,0)  in  the  rectangular  lattice.  Vein  dielectric  permittivity 
8.90 


structures.  Note  that  here  is  a  clear  phase  difference  between  the  first  two  rectangular  modes. 
Mode  plots  for  the  third  frequency  are  shown  in  Figure  25. 


For  the  third  modal  frequency,  the  square  lattice  is  dominated  by  an  array  of  dipoles 
oriented  in  the  x  direction.  The  magnitude  of  the  associated  fluctuations  is  quite  small  with 
respect  to  the  fundamental  mode.  The  dipoles  are  centered  within  the  air  gaps;  their  senses 
alternate  from  those  in  adjacent  cells.  On  the  other  hand,  E-  in  the  rectangular  lattice  is  composed 
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of  more  complicated  dipoles  oriented  in  the  y  direction.  There  is  also  a  periodicity  cast  in  the  x 
direction,  yet  it  extends  beyond  the  segment  of  the  lattice  shown.  The  peak  magnitudes  are 
centered  along  the  vertically  oriented  veins.  The  nodes  are  located  also  along  these  veins  but 
occur  at  the  air  gap  center.  The  peaks  are  juxtaposed  to  sides  of  the  center  making  for  an 
interesting  wave  form.  This  effect  is  a  sign  of  higher  order  modal  behavior.  Of  course,  there  are 
also  lines  of  vertically  oriented  nodes  as  indicated  at  the  right  side  of  the  rectangular  lattice  plot. 
This  behavior  is  contrasted  with  that  occurring  for  s  =  1 .52  in  the  next  section. 
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Figure  25.  Real  part  of  Ez  for  the  TM  mode  at  the  third  modal  frequency  computed  at  (left)  k-point  (wave  vector)  2 nla 
(0.5, 0.5,0)  in  the  square  lattice  and  at  (right)  k-point  2 nla  (0.5, 1.9,0)  in  the  rectangular  lattice.  Vein  dielectric  permittivity 
8.90 


4.4  Transverse  Magnetic  (TM)  Wave  Propagation  at  s  =  1 .52 

TM  electromagnetic  wave  propagation  has  been  analyzed  for  the  two  test  lattices  with  the 
index  of  refraction  set  to  1.52  in  the  dielectric  veins.  Band  diagrams  for  the  lattices  are  provided 
in  Figure  26.  In  this  case,  no  band  gaps  are  observed  in  either  lattice.  The  rectangular  lattice  does 
have  the  jagged  appearance  that  we  encountered  in  preceding  section.  We  continue  with  a  brief 
examination  of  the  electric  field  mode  shapes.  Figure  27  contains  plots  for  Ez  mode  shapes 
at  the  fundamental  frequencies  for  both  the  square  and  rectangular  lattices.  The  square  lattice  is 
dominated  by  monopoles  centered  at  the  vein  spots.  The  rectangular  lattice  is  characterized  by 
dipoles  oriented  in  the  y  direction.  The  peak  field  distributions  are  centered  on  the  horizontal 
dielectric  veins.  There  is  also  a  secondary  periodicity  in  the  x  direction  that  extends  beyond  the 
lattice  segment  shown.  The  field  magnitudes  differ  by  less  than  a  factor  of  two  between  the  two 
lattices  at  the  fundamental  frequencies.  In  Figure  28,  the  modal  shapes  for  Ez  are  plotted.  The 
square  lattice  is  dominated  by  a  complicated  mode  shape  with  nodes  distributed  along  the  veins. 
As  previously  encountered,  the  rectangular  lattice  is  dominated  by  dipoles  oriented  in  the  y 
direction,  but  with  a  long  wave  periodicity  in  the  x  direction.  In  this  case,  the  rectangular 
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Figure  26.  Dispersion  relations  for  transverse  magnetic  (TM)  wave  propagation  in  the  square  and  rectangular  dielectric 
vein  lattices  at  relative  permittivity  1.52.  The  dispersion  relation  for  the  square  lattice  is  shown  in  the  left  plot  while  the 
plot  for  rectangular  lattice  is  presented  on  the  right.  The  notation  wave  vector  number  is  used  to  indicate  the  location  of 
the  wave  vector  on  the  boundary  of  the  irreducible  Brillouin  zone.  Only  the  first  five  eigenmodes  are  displayed. 
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Figure  27.  Real  part  of  Ez  for  the  fundamental  TM  mode  computed  at  (left)  k-point  (wave  vector)  2 nla  (0.5, 0.5,0)  in  the 
square  lattice  and  at  (right)  k-point  2 nla  (0.5, 1.9,0)  in  the  rectangular  lattice.  Vein  dielectric  permittivity  1.52 

lattice’s  field  magnitude  exceeds  the  magnitude  on  the  square  lattice.  Figure  29  contains  the 
contours  for  the  real  part  of  E-  at  frequencies  for  the  third  mode.  This  electric  field  component 
again  has  a  very  compact  appearance  on  the  lattice.  Field  oscillations  are  centered  at  the  vein 
spots,  but  the  magnitudes,  at  least  with  respect  to  the  fundamental,  are  very  small.  The  sense  of 
the  field  reverses  at  alternating  spots.  In  some  sense,  the  field  behaves  in  a  dipolar  manner,  but 
the  fluctuations  are  complex.  On  the  other  hand,  the  rectangular  field  has  complicated  dipolar 
characteristics  with  a  primary  orientation  in  the  y  direction.  Peaks  appear  on  the  horizontal  veins 
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Figure  28.  Real  part  of  Ez  for  the  TM  mode  at  the  second  modal  frequency  computed  at  (left)  k-point  (wave  vector)  2n/a 
(0.5, 0.5,0)  in  the  square  lattice  and  at  (right)  k-point  2n/a  (0.5, 1.9,0)  in  the  rectangular  lattice.  Vein  dielectric  permittivity 
1.52 
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Figure  29.  Real  part  of  Ez  for  the  TM  mode  at  the  third  modal  frequency  computed  at  (left)  k-point  (wave  vector)  2n/a 
(0.5, 0.5,0)  in  the  square  lattice  and  at  (right)  k-point  2n/a  (0.5, 1.9,0)  in  the  rectangular  lattice.  Vein  dielectric  permittivity 
1.52 


and  at  every  third  of  the  air  gap.  Of  course,  there  is  a  long  wave  periodicity  that  extends  in  the  x 
direction  outside  of  the  lattice  segment  shown.  A  significant  investment  of  time  can  be  required 
to  fully  examine  and  understand  the  implications  of  both  TE  and  TM  modal  solutions.  It  is  also 
important  to  note  that  we  have  only  presented  the  shapes  for  three  distinct  frequencies  associated 
with  each  of  these  modes.  Obviously,  there  are  many  other  modes.  Instead  of  exploring  higher 
order  modes,  it  is  instructive  to  turn  our  attention  to  interpreting  the  impact  on  light  coloration 
associated  with  some  of  the  modes  detailed  above. 
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4.5  Comparing  Results  for  the  Two  Lattices 

In  preceding  section,  we  have  discussed  a  modal  analysis  of  electromagnetic  wave 
propagation  for  square  and  rectangular  2D  crystal  lattices.  Our  interest  now  lies  in  studying  how 
differently  light  propagates  on  the  rectangular  lattice  than  on  the  square  lattice.  For  high  relative 
dielectric  permittivity  values,  the  square  lattice  possesses  band  gaps.  Since  the  rectangular  unit 
cell  geometry  is  “closer”  in  shape  to  the  photonic  crystal  mentioned  in  Section  3,  one  may  ask  if 
band  gaps  can  appear  for  this  lattice,  especially  for  the  relative  dielectric  permittivity  of  chitin,  a 
value  quite  near  unity.  The  square  lattice’s  field  exhibits  highly  uniform  mode  shapes  when 
sampled  at  the  so-called  M  k-point  2n la  (0.5. 0.5.0).  The  electric  field  fluctuations  are  usually 
centered  on  each  lattice  spot,  air  gap  center  or  at  the  midpoint  of  the  dielectric  vein.  There  are 
phase  differences  between  the  x  and  y  components  of  the  field  lattice,  but  the  phase  difference  is 
the  width  of  one  unit  cell.  This  behavior  is  echoed  for  the  square  lattice  at  the  lower  value  of 
s  =  1 .52 ,  yet  in  this  case,  no  band  gaps  exist.  Band  gaps  are  usually  created  when  the 
differences  in  s  are  high  across  the  lattice  since  the  field  concentrates  in  regions  of  high  s  .[2] 
At  this  reduced  dielectric  value,  the  differences  are  too  small  to  generate  band  gaps. 

Our  interest  is  to  see  if  similar  behavior  occurs  for  the  rectangular  lattice,  even  at  the 
higher  dielectric  value.  The  calculations  show  that  the  rectangular  field  is  very  different. 
Although  the  field  fluctuations  are  periodic  on  the  rectangular  lattice,  the  periodicity  in  the  x 
direction  is  manifested  over  a  distance  of  several  unit  cells.  That  is  to  say,  for  the  square  lattice, 
the  period  for  a  complete  field  oscillation  is  no  more  than  two  cells.  For  the  rectangular  case,  the 
period  in  the  y  direction  is  usually  two  cells,  but  in  the  x  direction,  a  period  of  field  oscillation 
may  encompass  eight  unit  cells.  Also,  on  the  rectangular  lattice,  the  field  components  usually 
consist  of  dipoles  that  are  elongated  in  the  x  direction.  As  the  modal  order  increases,  these 
dipoles  tend  to  adopt  more  complex  structure.  On  the  other  hand,  the  square  lattice’s  field  adopts 
complexity  through  both  dipole  and  quadrupole  behavior  centered  in  the  each  cell,  vein  spot  or 
vein  midpoint.  There  is  also  a  phase  difference  observed  for  the  rectangular  lattice,  but  in  the  x 
direction,  the  difference  is  comprised  of  several  cell  widths.  The  reason  for  the  field  behavior 
observed  for  the  rectangular  lattice  is  likely  due  to  field  concentration  caused  by  the  close 
spacing  of  the  dielectric  veins  in  the  y  direction.  Recall  that  a  unit  cell  is  about  26%  as  wide  as  it 
is  long.  Since  the  electric  field  tends  to  reside  mostly  within  the  higher  dielectric  material,  the 
close  spacing  tends  to  elongate  the  field  in  the  x  direction  making  band  gap  formation  difficult, 
especially  for  lower  values  of  s .  It  is  necessary  to  get  the  dielectric  and  air  bands  to  separate  in 
order  for  band  gaps  to  form.  It  simply  does  not  happen  for  the  rectangular  lattice  under  these 
conditions.  Still,  it  is  important  to  say  that  this  assertion  is  merely  conjecture,  not  proof. 

4.6  Example  Calculation  Set  at  the  Micron  Scale 

Preceding  discussions  have  presented  a  modal  analysis  of  electromagnetic  wave 
propagation  within  dimensionless  2D  square  and  rectangular  vein  crystal  lattices  for  two 
different  relative  dielectric  permittivity  values.  In  this  section,  we  apply  dimensions  to  the 
rectangular  model  and  deduce  the  wavelengths  of  light  propagated  for  different  wave  vectors 
existing  at  the  border  of  the  Brillouin  zone.  The  selection  of  a  particular  k-point  implies  a 
propagation  direction  for  the  mode.  Since  mirror  symmetry  is  enforced  in  this  study,  all  modes 
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(electromagnetic  waves)  are  confined  to  the  xv-plane.  The  wavenumber  k  for  the  mode  is  simply 
derived  as  the  magnitude  of  the  wave  vector,  i.e., 

k=  I  k  1=  Jk[+kf+kf  (49) 

From  the  wavenumber,  we  can  easily  determine  the  wavelength  A  for  the  mode.  To  do  so,  we 
use  the  formula 


(50) 


For  dimensionless  input,  the  MPB  computer  code  produces  dimensionless  output.  [2]  This 
procedure  offers  great  utility  since  dimensionless  solutions  can  be  scaled  into  any  system  of  units 
desired.  The  dimensionless  wave  vector  k  may  be  expressed  in  terms  of  the  dimensioned  wave 
vector  k  as  follows 


k  = 


k  a 
2 ~k 


(51) 


The  dimensionless  angular  frequency  co  is  expressed  in  terms  of  the  dimensioned  angular 
frequency  by  the  formula 


_  coa 

O)  = - 

2.71  c 


(52) 


where  c  is  the  speed  of  light  in  vacuum,  approximately  2.997  x  1 08  m/s.  Equations  (51)  and  (52) 
may  be  used  to  convert  k  and  co  into  their  dimensioned  equivalents. 

For  a  specific  example  regarding  the  lattice  shown  in  Figure  6,  the  lattice  constant  is  set 

to  0.333  pm  (or  0.333  xl  0”6  m).  By  sweeping  over  the  wave  vectors  at  the  edge  of  the  Brillouin 
zone  shown  in  Figure  7,  we  obtain  the  range  of  wave  numbers  and  wavelengths  shown  in  Figure 
30  for  s  =  8.90 .  As  is  evidenced  by  Figure  30,  for  the  selection  of  wave  vectors,  the  rectangular 
lattice  mostly  admits  light  waves  in  the  range  of  0.25  pm  to  3.5  pm  with  corresponding  wave 
numbers  between  0  rad/pm  to  40  rad/pm.  The  wavelength  range  is  of  course  broader  than  that 
stated  due  to  the  singularity  at  k  =  0.  Recall  that  the  Brillouin  zone  does  include  the  wave  vector 
k  =  (0,0,0) .  The  frequency  loci  can  also  be  transformed  into  a  common  system  of  units  for  each 
mode.  Figure  31  contains  this  plot  for  s  =  8.90 . 
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Figure  30.  Wavelengths  and  wave  numbers  predicted  for  a  rectangular  lattice  with  lattice  constant  set  at  0.333  pm  and  a 
vertical  aspect  ratio  of  0.26315.  The  wavelength  unit  is  micrometers;  the  wave  number  unit  is  tens  of  radians  per 
micrometer.  Relative  dielectric  permittivity  is  set  at  8.90 
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Figure  31.  Band  or  mode  frequencies  predicted  for  a  rectangular  lattice  with  lattice  constant  set  at  0.333  pm  and  a 
vertical  aspect  ratio  of  0.26315,  relative  dielectric  permittivity  8.90.  The  frequency  unit  is  Petahertz. 

The  dimensioned  frequency  plot  has  an  appearance  that  is  identical  to  that  for  the  dimensionless 
frequency  graph  shown  in  Figure  8.  Now  real  frequency  units  are  assigned.  Similar  calculations 
can  be  performed  on  the  same  lattice  with  s  =  1 .52 .  The  resulting  wavelength  and  wave  number 
distributions  are  the  same  as  shown  in  shown  in  Figure  30  since  the  Brillouin  zone  is  the  same. 
The  attendant  frequency  plot  is  provided  in  Figure  32.  On  the  whole,  a  major  effect  of  reducing 
the  relative  permittivity  on  the  rectangular  lattice  is  to  elevate  the  modes  to  higher  frequencies 
and  to  increase  the  jaggedness  of  the  frequency  plots.  In  some  cases,  the  separation  in  frequency 
between  modes  also  shrinks.  In  terms  of  light  coloration,  this  example  largely  involves  light  in 
the  wavelength  range  of  0.25  pm  to  3.5  pm  (250  nm  to  3,500  nm).  For  the  first  Brillouin  zone, 
this  band  is  roughly  concurrent  with  the  visible  spectrum  and  the  near  infrared.  As  one  may 
conclude,  this  solution  procedure  renders  a  wealth  of  information  on  the  behavior  of 
electromagnetic  radiation  in  the  photonic  crystal  and  allows  very  precise  control  over  wave 
number  and  frequency  combinations  for  research  consideration. 
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Figure  32.  Band  or  mode  frequencies  predicted  for  a  rectangular  lattice  with  lattice  constant  set  at  0.333  pm  and 
vertical  aspect  ratio  of  0.26315,  relative  dielectric  permittivity  1.52.  The  frequency  unit  is  Petahertz. 
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5  CONCLUSIONS 

A  frequency-based  solution  method  for  the  propagation  of  electromagnetic  waves  has 
been  discussed  in  this  report.  Specifically,  modal  frequencies  and  shapes  for  square  and 
rectangular  2D  photonic  crystals  have  been  calculated  for  two  different  relative  dielectric 
permittivity  values.  In  this  study,  the  high  permittivity  values  are  confined  to  dielectric  veins.  A 
principal  result  is  that  for  the  high  permittivity  case,  a  band  gap  exists  for  the  square  lattice  in 
both  the  transverse  electric  and  magnetic  waves.  No  such  gaps  exist  for  the  specific  rectangular 
lattice  (rectangular  ratio  -  1:0.26)  considered  here.  For  the  low  permittivity  analyses,  no  band 
gaps  exist  for  either  lattice.  These  calculations  may  be  repeated  for  the  higher  permittivity  value 
s  =  2.31  in  direct  correspondence  with  a  proper  value  for  chitin.  However,  this  oversight  has 
little  impact  on  the  results.  The  frequencies  are  affected  only  to  a  small  degree  for  the  higher 
order  modes.  Details  concerning  the  electric  field  components  have  been  discussed  for  each 
lattice.  Significant  differences  between  electric  fields  for  the  square  and  rectangular  lattices  have 
been  highlighted.  The  level  of  field  compression  caused  by  the  close  spacing  of  dielectric  veins 
seems  to  cause  waves  on  the  rectangular  lattice  to  depart  from  those  observed  for  the  square 
lattice.  A  point  of  future  interest  may  be  to  determine  what  minimal  changes  need  be  made  to  the 
rectangular  lattice  in  order  to  motivate  the  appearance  of  band  gaps.  We  may  also  continue  to 
add  geometric  complexity  to  the  rectangular  lattice  in  order  to  make  it  better  representative  of  the 
butterfly  wing  lattice  introduced  early  in  this  report.  Such  investigations  are  under  consideration 
for  a  forthcoming  program  of  research.  As  a  collateral  benefit,  this  project  has  served  to  exercise 
the  author’s  skills  both  in  utilizing  MPB  in  working  with  the  equations  of  electromagnetism. 
Musculature  for  the  latter  has  sadly  atrophied  a  bit  with  the  passage  of  time. 
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